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Abstract 

Bursty transport phenomena associated with convective motion present universal statistical char- 
acteristics among different physical systems. In this letter, a stochastic univariate model and the 
associated probability distribution function for the description of bursty transport in plasma tur- 
bulence is presented. The proposed stochastic process recovers the universal distribution of density 
fluctuations observed in plasma edge of several magnetic confinement devices and the remarkable 
scaling between their skewness S and kurtosis K. Similar statistical characteristics of variabili- 
ties have been also observed in other physical systems that are characterized by convection such 
as the X-ray fluctuations emitted by the Cygnus X-l accretion disc plasmas and the sea surface 
temperature fluctuations. 
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Plasma turbulence and the associated heat and particle transport play a major role in the 
levels of plasma confinement in magnetic fusion devices. In plasma edge the turbulent fluc- 
tuations are large and bursty due to the convective motion of strongly nonlinear structures 
formed during the nonlinear saturation of plasma instabilities. Such coherent structures in 
an unambiguous manner effectively contribute to radial transport and to intermittency [l|]. 
As a result, the transport process departs from the diffusive picture associated with the 
Gaussian case of weak independent fluctuations. 

In order to understand the underlying mechanism of the turbulent transport in the plasma 
edge, it is crucial to investigate the statistical characteristics of these bursty fluctuations. Ex- 
perimental investigations have indeed revealed the bursty nature of particle transport in the 
scrape-of-layer (SOL) of magnetically confined plasmas [2j. The appearance of structures 
such as plasma blobs, avaloids [3| is attributed to the formation of field-aligned structures - 
induced by the charge separation of the magnetic curvature drifts - that propagate radially 
far into the SOL. 

The statistical behavior of the density fluctuations associated with such bursty dynamics 
has been investigated in several experiments. A comprehensive study |4j , that included mea- 
surements from various magnetic confinement devices - including TORE SUPRA, MAST 
and ALCATOR C-MOD - showed that the associated extreme probability distribution func- 
tions (PDFs) are universal (see Fig. [1]) in the sense that have the same properties in many 
confinement devices with different configurations. The investigations were carried out using 
Langmuir probes that measure the saturation current. 

Remarkably, similar [5[| extreme distributions have been also observed for the bursty X-ray 
fluctuations associated with transport events in the Cygnus X-l accretion disc plasmas that 
are linked to instabilities which give rise to turbulent transport and extreme statistics (see 
Ref. [6] and Fig. 12 therein). However, known forms of extreme PDFs - such as the Frechet 
or the Gumbel distributions - do not have [s] the proper form to fit well the experimentally 
observed distributions associated with bursty convective transport processes. 

Bursty dynamics is characterized by strongly non-Gaussian PDFs and non-vanishing 
probabilities of extreme events. In such dynamical systems, the higher order moments are 
commonly used to determine the scaling properties of the fluctuating fields. The non- 
Gaussian features are usually quantified in terms of skewness and kurtosis of the PDF of 
fluctuating fields. For a centered random variable x, i.e. (x)=0, with variance a 2 , skewness is 
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FIG. 1: The PDF plot of the ion saturation current in the Tore Supra (solid line), Alcator C-Mod 
(thick solid line), MAST (dashed-dotted line), and PISCES (dots). The ion saturation current was 
normalized to the standard deviation and the integral of the four PDFs is set equal to 1. Figure 
reprinted with permission from Ref. J] (http://link.aip.Org/link/7PHPAEN/10/419/l), Fig. 3, 
Copyright 2003, American Institute of Physics. 



defined by the S = (x 3 ) /a 3 and the kurtosis (in the newer literature) by K = (x 4 ) /a A — 3. 
Skewness is a measure of asymmetry of a PDF; if the left tail is more pronounced than 
the right tail, the PDF has negative skewness and when the reverse is true, it has positive 
skewness. Kurtosis measures the excess probability (flatness) in the tails, where excess is 
defined in relation to a Gaussian distribution. For Gaussian distributions both S and K are 
equal to zero. 

Using ten thousand observed density fluctuation signals measured in TORPEX, Labit et 
al. [t| showed that a unique parabolic scaling relation holds, K = 1.5025* 2 — 0.226, between 
the skewness S and the kurtosis K (see Figure 1 of Reference 7|). It was also shown that 
the PDFs of the measured signals, including those characterized by a negative skewness 
can be described by a special case of the Beta distribution. The density fluctuations were 
associated with regimes of drift-interchange (D-I) turbulence generated in regions of bad 
magnetic field curvature and convected away by the E x B fluid motion. 

Remarkably, there is a striking similarity of the observed K—S scaling with that of sea 
surface temperature (SST) fluctuations that are governed by advection through ocean cur- 
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rents [9] . Sura and Sardeshmukh j^] proposed a nonlinear Langevin model that can predict 
the observed scaling in some limits of the model parameters while Krommes [8| generalized 
it to include linear waves - an essential feature of D-I turbulence - and he numerically cal- 
culated the associated PDF that includes four indepedent parameters. The fact that the 
same scaling applies to different physical situations leads to conjecturing that the scaling 
arises due to basic constraints. The general mathematical constraints on the K — S relation, 
as arise from the definition of kurtosis and skewness jlO 1 ], do not provide any insight about 
the observed scaling. A recent comprehensive study on the observed K — S scaling between 
various fusion devices showed that the data align along parabolic curves 11]. A phenomeno- 
logical model using the assumption that the fluctuating signals include a linear combination 
of two basis (Gamma or Beta ditributions) PDFs attempting to accommodate experimental 
evidence was also provided. However, the large number of the free parameters entering to 
the model and the small difference between using one or a sum of two Beta PDFs make the 



model less flexible 
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It becomes evident that these observations among different physical systems reveal a uni- 
versal character associated with strongly non-Gaussian processes. However, the key question 
is what kind of underlying mechanism is responsible for the universally observed statisti- 
cal features and what PDF can describe them. In this Letter, a univariate model for the 
statistical description of bursty fluctuations is presented, using as an example the aforemen- 
tioned plasma edge density fluctuations. Furthermore, the associated PDF is derived and 
it is shown that it recovers the universally observed distributions, documented in Ref. 4j 
and the remarkable parabolic scaling between skewness and kurtosis as well, documented 
in Ref. Q. The derived results have universal character, and thus may be applicable to all 
aforementioned physical systems. 

The non-linear processes described by the standard models of turbulence in magnetized 
plasma are quadratic and are linked with small or large scale convection processes attributed 
to electric drifts. Thus, it is natural one to assume that the universally observed statistical 
characteristics associated with the bursty behavior of fluctuations may be attributed to 
processes that emerge from the non-linear quadratic interaction between turbulent fields. 
However, extreme statistical features are expected to appear when strong non-Gaussian 
processes coexist with Gaussian ones. In order to describe the associated universal statistical 
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FIG. 2: (a) Example of bursty time-series W(t) for different values of 7: (a) -0.11, (b) 0.031, (c) 
0.98, corresponding to negative (top), positive (middle), and significantly large positive (bottom) 
skewed distributions respectively. 



properties, we propose a univariate non-Gaussian process W(t) given by: 

Z W -<Z( t) ) Z»(«) - <*(«)) 

which results from the superposition of a Gaussian Z(t) (with standard deviation er ) with 
a non-Gaussian process attributed to the square of Z(t). The latter corresponds to the 
strongest quadratic non-linearity that may arise due to the interaction of fluctuating fields. 
The coefficient 7, in front of the quadratic non-Gaussian component is a parameter that 
measures the deviation of W(t) from Gaussianity. The process W(t) may well describe 
extreme bursty behavior of fluctuations that is characterized by non-linear structures (blobs, 
avaloids) that travel (convect) through a sea of Gaussian fluctuations. Figure [2] presents 
some typical time series W(t) calculated from Eq. (1) by using a centered random Gaussian 
fluctuation Z(t). In all cases the bursty nature of W(t) is evident. For the sake of simplicity, 
from now on, we drop the dependence of Z and W on time and consider that Z is a centered 
Gaussian process. It should be noted here that the process described by Eq. (1) has been 
presented in the literature 12| as a characteristic simple-to-construct non-Gaussian process 
and was used as an example for the calculation of the systematic errors of covariances and 
moments up to fourth order of non-Gaussian time series. 

For the derivation of the PDF p(w) of the random process W(t), we express the cumulative 
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distribution function (CDF) P(w) as follows 

P(w) = Pr (W < w) = Pr (7 (Z - Zl ) (Z - z 2 ) < 0) 



(2) 



where z\ and z 2 are the roots of the polynomial 



/(Z)=7^=^ + --™. (3) 
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Equation [2] can be expressed in terms of the CDF P(Z) at z = Z\, z 2 and p(w) can be 
derived by a simple differentiation that leads to the following expression: 

p ( w j _ J °' when ^ W = 47( ^ W - ^o) < else; 

P " i cosh (^) exp (-^) ' 

where wq = — (47) _1 —7. The value and the sign of parameter 7 controls the shape and the 
range of non-zero values of p(w). For 7 = 0, the p(w) reduces to that of a centered Gaussian 
random process. For 7 > (7 < 0) a cut off value wo exists and p(w) gets asymmetric 
presenting long tails in the positive (negative) axis (cf. Figure 3). The minimum (in absolute 
sense) value of the cut-off is equal to w 0m = ±1 and corresponds to 7 = =Fl/2. Note, that 
the same analytic expression for p(w) can also be derived by noticing that the non-Gaussian 
process W can be re-written in terms of a scaled and shifted non-central chi-squared random 
process with one degree of freedom, 

( Z 1 \ 2 

W = w + j[ — + — ) (5) 
V^o 277 

and using the expression of non-central chi-squared PDF along with standard transform 
techniques of random variables. 

In Fig. HI we have plotted p(w) choosing positive values of 7 that range up to 7 = 0.2. 
The resulting distributions exhibit the same behavior with the measured distributions of 
plasma edge density fluctuations, as have been observed in several magnetic confinement 
devices and presented here in Fig. [TJ Note, that the characteristic clustering of the PDF 
curves around 0.05 for values around 2.5 which is experimentally observed (cf. Fig. 1) is 
recovered by the distribution p(w) (cf. Figs. 3 and 4) - independently on the value of 7 > 
- a characteristic that is not recovered in the plots of PDF in Ref. Similar features 
appear also in the PDF of the x-ray observations associated with anomalous transport in 
accretion disks (cf. Fig. 8 of Ref. |6j). 
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FIG. 3: The family of distributions p(w), for different values of 7 = —2, —0.02, 0, 0.02 and 2. The 
symmetry with respect to 7 is evident. As 7 increases the absolute value of the cut-off decreases 
(increases) for 7 < 0.5 (7 > 0.5). 
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FIG. 4: The family of distributions p(w), for 7 = 0.05, 0.1 and 0.2. The resulting PDFs recover 
the universally observed distributions presented at Fig. [TJ 

Unlike, to the BHP distribution [13] which does not have any free parameter, p(w) de- 
pends on 7 such that its higher order moments can receive multiple values. Furthermore, the 
existence of a single free parameters allows p(w) to be used for the fitting of experimentally 
observed distributions. The high order moments of multivariate Gaussian processes can be 
determined by a simple method based on the Wicks theorem. All odd moments are zero and 
all even moments can be reduced to homogeneous polynomials. For the considered process 
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W, the values of skewness and kurtosis depend on the value of 7 and are given by [121 ]. 

3 + 4 7 2 4 !+T 2 

For 7 = 0, the process W is Gaussian and K = S = 0. As | 7 | increases kurtosis and skewness 
converge to their extreme values K max = 12 and S max = ±2 3 / 2 , respectively. Noticeably, 
these values agree with the range of the values reported in the Letter by Garcia et al 
[yj] , in which the authors investigate intermittent transport in plasma edge using extensive 
numerical simulations that take into account the presence of SOL. 

Using the parabolic ansatz, K = aS 2 + b for the parametric relations in Eq. (jSJ), it can 
be easily found that 

nM 1 o (l + T 2 )(l + 27 2 ) m 
° (7) = 12 (3 + 4 7 ^ ' (?) 

and 6 = 0. The function 0(7) converges rapidly to the value a = 3/2,which is exactly the 
value of the universally observed parabolic scaling, reported in Refs. [3, la, |9j . 

It is interesting to compare the derived K — S curve of Eq. ([6]), to that attributed to a 
quadratic product of two central Gaussian processes Zi, Z 2 , i.e. Wr = r YrZiZ 2 . Here 7 r 
denotes the correlation between the Gaussian fields. For Z% = Z2 and 7 r = 1, the PDF 
of Wr is simply the chi-squared distribution with one degree of freedom. The associated 
PDF p(Wr) has been presented in Ref. [151 ]. and the corresponding values of skewness and 
kurtosis were found equal to: 

q 9 3 + 7r , „ , l + 6 7 g + 7^ ,^ 
Sr = - 27r (l + 7 g)3/ 2 and R r = 6 (i + 7 g)2 . ( 8 ) 

respectively. The parabolic relation, K = arS 2 + b r , between Sr and K T) results to the 
derivation of the following coefficients: 

fl r = and b T = 6. (9) 

Note that for | 7 r| = 1, the values of kurtosis and skewness are equal to K max and S max , 
while for 7r = 0, are equal to 6 and 0, respectively. 

Equations and (jUJ) define a closed curve in the K — S space for continuous values of the 
parameters 7 and 7r respectively (see Fig. [5]). The low boundary of the curve corresponds 
to extreme bursty processes W given by Eq. (1) and distributed according to p(w), while 
the upper boundary corresponds to Wr which is distributed according to the "local flux 
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FIG. 5: The K — S curve for stochastic processes described by quadratic polynomial of gaussians. 
The solid line corresponds to W(t) and the dashed line to wr(t) process. 



PDF" (see Eq. (7) in Ref. |15|). The clustering of the K — S values around the K = 1.5S 2 
curve is patently seen in Fig. 1 of Reference 7(] and in Fig. 3 of Reference 9] . 

It is straightforward to show numerically, that for non-Gaussian processes (described by 
second order polynomials of Gaussians) the "cloud" of K — S points fall into (not shown 
here) the closed curve. The latter follow a parabolic trend - similar to the experiments 11] 
- with scaling that depends on the selection of the polynomial coefficients. Furthermore, for 
processes described by higher order Gaussian polynomials the "cloud" follows also parabolic 
trend including much higher values of K and S. 

In conclusion, a univariate stochastic model for the description of extreme bursty fluctu- 
ations W based on generic properties of quadratic non-linearities is presented. The model 
and the associated probability distribution function describe statistical properties of the SOL 
density fluctuations which are governed by the E x B convection. The universal character of 
the stochastic process stems from the fact that the associated extreme PDF recovers proper- 
ties of variabilities (density/sea temperature/x-ray intensity) observed at boundary regions 
(SOL/sea surface/accretion edge) of different physical processes that are characterized by 
convection (electric drift/ocean current/rotation). The proposed univariate stochastic model 
describes the statistical characteristics of these relaxation phenomena at a state of extreme 
statistical behavior. It is evident that the parabolic relation between S and K when observed 
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provide relevant information about the underlying processes. 
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